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A  Cluster  Analysis  Program  For  Image  Segmentation 
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M.  F.  Janowitz 


I.  Introduction. 


A  more  detailed  description  of  the  contents  of  this 
paper  occurs  in  and  £8] .  The  present  version  is  intended 
only  to  be  a  survey,  and  for  that  reason  will  omit  proofs  of 
results.  The  thrust  of  the  work  is  to  relate  the  image 
segmentation  problem  to  a  rather  general  model  for  cluster 
analysis  that  was  introduced  in  C3f]»  Some  techniques 
suggested  by  the  model  are  described,  and  their  implementation 
illustrated.  Section  2  is  devoted  to  some  background 
material  from  the  theory  of  partially  ordered  sets,  with 
Section  3  containing  the  image  segmentation  model.  Section  4 
has  an  elementary  discussion  of  some  of  the  underlying 
statistical  considerations,  with  Sections  5  and  6  devoted  to 
the  description  of  certain  segmentation  techniques  suggested 
by  the  model.  Finally,  in  Section  7  the  techniques  are 
implemented  on  real  data,  and  the  remits  are  presented. 
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Background  material  from  the  theory  of  partially  ordered 


sets. 

A  basic  familiarity  with  the  theory  of  partially  ordered 
sets  and  lattices  will  be  assumed.  Despite  this,  it  will  be 
convenient  to  specifically  develop  certain  terminology 
here.  Unless  otherwise  specified,  all  partial  orders  will 
be  denoted  4:  ,  with  reserved  for  set  inclusion.  The 
symbol  <p(X)  represents  the  Boolean  algebra  of  all  subsets 
of  the  set  X,  ordered  by  set  inclusion.  Unions  and  inter¬ 
sections  have  the  standard  symbols  of  0  and  O  ,  and  it 
will  be  useful  to  let  R  denote  the  nonnegative  real  numbers, 
ordered  in  the  usual  manner. 

Let  P,  Q  be  partially  ordered  sets.  A  mapping  ^  :P-*Q 
is  called  iso  tone  if  p  ^  p2  in  P  implies  that  $(p  )  -  4>(P2) 

in  Q.  It  is  residuated  if  it  is  isotone  and  there  i3  an 
isotone  mapping  -*  P  such  that  p  4s  <$Hl>(p)  and  q  ^  $>4>+(q) 

for  every  p  e  P,  q  e  Q.  The  mapping  b+  is  the  residual 

mapping  associated  with  <p  ,  and  it  is  completely  determined 

by  4);  likewise,  <p  is  completely  det  mined  when  4>+  is 

specified,  so  for  a  given  residual  mapping  0  ,  it  will  often 

*  .... 
be  convenient  to  write  0  for  the  residuated  mapping  with  which 

0  is  associated.  The  notation  Res(P,Q)  or  Res(P)  when  P  =  Q 

will  be  used  to  denote  the  set  of  all  residuated  mappings 

from  P  into  Q,  with  Res+(Q,P)  and  Res+(P)  for  the  corresponding 

sets  of  residual  mappings.  It  turns  out  that  residuated 

and  residual  mappings  play  a  natural  role  in  the  Jardine- 

Sibson  model  C9J  for  hierarchical  clustering.  This  is 

because  ( C3]  ,  Lemma  4.1,  p.60)  C:  R-»C?(X) ,  where  X  is  a 

finite  nonempty  set,  is  residual  if  and  only  if  it  satisfies 

the  following  three  conditions: 
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(1)  C  is  isotone. 

(2)  C(h)  =  X  for  some  h  e  R. 

(3)  Corresponding  to  each  h  £  R.  there  is  a  positive 
real  number  6  =  6(h)  such  that  C(h)  =  C(h  +  6). 

If  X  is  the  set  of  all  2  element  subsets  of  the  set  P, 
this  turns  out  to  be  equivalent  to  C  being  a  numerically 

stratified  clustering  in  the  sense  of  Jardine  and  Sibson 

00  /  P-61. 

The  theory  of  residuated  mappings  is  rather  extensively 
developed  in  OQ  ,  and  the  reader  is  referred  to  that  source 
for  further  information. 


Ill .  The  image  segmentation  problem. 

The  underlying  input  data  is  a  function  F:  P-*R, 
where  P  represents  a  closed  bounded  rectangle.  One  thinks 
of  P  as  representing  a  picture  having  certain  natural  regions 
that  are  characterized  in  some  known  manner  by  the  values 
of  F.  This  might  involve  regions  on  which  F  is  constant,  or 
some  sort  of  texture  measure,  or  some  combination  of  attrib¬ 
utes  of  F.  Unfortunately,  no  direct  knowledge  of  the  values 
of  F  is  available.  Rather,  the  actual  input  is  a  function 
G:  X-*  R,  where  X  =  {l,2,...m}  X  {l,2,...,n}  for  suitable 
positive  integers  m  and  n.  One  can  think  of  G  as  being  some 
sort  of  sample  of  F  or  as  somehow  summarizing  the  values  of 

F,  possibly  by  averaging  F  over  some  small  subregions  of  P. 
Often  there  is  some  noise  involved  in  the  passage  from  F  to 

G,  so  G  can  only  be  regarded  as  an  estimate  of  F.  The  idea 
is  to  use  G  to  somehow  recapture  the  principal  regions  of 

P  as  they  are  defined  by  F.  A  more  precise  and  detailed 


description  of  these  underlying  assumptions  occurs  in  D1  , 
but  the  present  version  is  sufficient  for  the  current  dis¬ 
cussion. 

A  rather  broad  model  for  hierarchical  clustering  was 
presented  in  [3] .  It  included  both  the  Jardine-Sibson  model 
and  a  graph-theoretic  model  due  to  Matula  0-0]  .  It  will 
now  be  shown  that  it  also  includes  the  image  segmentation 
problem  as  we  have  stated  it.  Recall  that  the  actual 
input  data  is  a  mapping  G:  X“*R,  where  X  =  {1,2,..., m}  X 
{l,2,...,n}.  The  idea  is  to  transform  G  into  a  function 
H:  X->R,  where  H  has  fewer  distinct  levels  of  range  values 
than  does  G.  The  values  of  H  should  be  thought  of  as  an 
estimate  of  the  regions  represented  by  F.  There  is  a 
natural  bijection  between  mappings  G:  X-*R  and  elements  of 
Res'*  (R,C?(x) }  .  It  is  given  by  associating  with  G  the  unique 
residuated  mapping  G  :  (?(X)-»R  to  which  it  extends,  the 

*  1 

extension  being  given  by  the  rule  G  (M)  =  V[G(m)  :  m  t  . 

* 

One  then  takes  the  residual  mapping  G  associated  with  G  . 
Thus  the  image  segmentation  problem  may  be  viewed  as  the 
study  of  transformations  of  Res  +  (R,  (?(X) )  into  Res+(R,<?(Y) )  , 
where  Y  S  x.  Such  transformations  v  il  be  called 
segmentation  methods.  The  reason  fcr  Y  £  X  will  be  apparent 

from  the  concrete  examples  that  will  be  presented.  This 
places  the  segmentation  problem  squarely  within  the  frame¬ 
work  of  the  model  contained  in  C3). 

No  assumption  has  thus  far  been  made  as  to  whether  the 
input  data  has  ordinal  or  numerical  significance.  Suppose 
that  the  input  data  has  only  ordinal  significance  in  that 
the  actual  numerical  values  of  G  have  no  significance, 
but  one  might  still  be  able  to  attach  some  significance  to 
an  assertion  of  the  form  G(i,j)  <  G(k,l).  Is  it  possible 
to  describe  the  class  of  segmentation  methods  that  are  suit¬ 
able  for  this  type  of  data?  Indeed  it  is,  but  before  the 


description  can  be  presented,  some  additional  terminology 

is  required.  Suppose  one  is  given  a  segmentation  method 

and  a  residual  mapping  9  on  R.  To  say  that  ^3- is  9  — 

compa table  (  [3],  p.68)  is  to  say  that  for  every  C  e  Res* (RfG^X) ) , 

it  is  true  that  9-(C®0)  =  3-(C)®  6.  A  minimal  requirement  for 
3-  is  that  it  be  9-compatable  for  all  order  automorphisms 
6  of  R.  Such  techniques  are  called  monotone  equivariant  by 

Jardine  and  Sibson.  A  precise  mathematical  proof  is  given 
in  [53  of  the  fact  that  every  segmentation  technique  suitable 
for  use  with  ordinal  data  is  in  a  sense  equivalent  to  a 
monotone  equivariant  technique.  The  equivalence  is  the 
fact  that  if  the  outputs  were  rank  ordered,  they  would 
be  identical.  This  now  makes  one  inquire  into  the  nature 
of  monotone  equivariant  techniques.  This  question  was 
posed  and  answered  in  [43  (Theorem  1,  p.149) .  For 

C  e  Res  (R,(?(X) )  ,  one  says  that  h  (h  e  R)  is  a  splitting  level 

* 

of  C  in  case  h  is  the  image  of  some  subset  of  X  under  C  , 
the  residuated  mapping  with  which  C  is  associated.  If 
0  <  h1  <  ...  <  ht  denotes  the  sequence  of  splitting  levels 

of  C,  then  the  segmentation  method  3-  is  monotone  equivariant 
if  and  only  if  it  is  true  that: 

(i)  every  splitting  level  of  3(C)  is  a  splitting  level 
of  C, 

(ii)  the  sequence  9C(0)  '<  ScUi^)  <  . ..  £  3C(hfc)  =  Y 
depends  only  upon  the  sequence  C(0)  <  C(h^)  <  ... 

<  C (hfc)  =  x  and  is  independent  of  the  actual  values 
of  the  h^'s; 


(iii)  conditions  (i)  and  (ii)  hold  for  every  C  in 
Res*(R ,<?(X) )  . 


One  can  associate  with  each  subset  T  of  Res"MR)  the 
collection  a(T)  of  all  segmentation  methods  that  are  9- 
compatable  for  every  0  e  T.  If  the  nature  of  the  input 
data  makes  compatability  with  a  certain  set  T  of  residual 
mappings  desirable,  then  the  determination  of  a(T)  will 
produce  a  description  of  the  most  general  type  of  segment¬ 
ation  technique  that  one  should  properly  use.  We  have  just 
noted  that  if  T  denotes  the  set  of  all  order  automorphisms 
of  R,  then  a(T)  is  the  set  of  all  monotone  equi variant 
mappings .  An  extremely  important  result  is  provided  oy 
F.  Baulieu  [1]  .  It  yields  an  explicit  description  of  all 

classes  a(T)  that  can  arise  for  T  containing  the  set  of  all 
order  automorphisms  of  R.  Here  is  a  description  of  these 
classes: 


(51)  Flat  methods.  *JC(h)  depends  only  on  C(h). 

(52)  Semiflat  methods.  £C(h)  depends  upon  C(h)  and  C(0) 

(53)  Divisive  methods.  ^Cfh^)  depends  upon  C(h^), 
C(h^+1),  ...»  C(ht),  whe’"'  hfc  is  the  highest 
splitting  level  of  C. 

(54)  SF  (1)  methods.  fJC(O)  depends  only  on  C(0); 

for  h^  0,  C(h^)  depends  upon  C(h^) ,  C(hi+1), 

•  •  •  t  C(n^). 

(55)  Monotone  equi variant  methods. 

The  simplest  of  these  classes  of  segmentation  methods 
is  of  course  the  class  of  flat  methods,  and  it  is  to  this 
class  that  we  now  direct  our  attention.  Such  methods  are 
easy  to  describe  and  easy  to  implement.  Given  a  fixed  non- 
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empty  finite  set  X  and  an  increasing  sequence 
MiC  M2C  •••  C  Mt  =  X 

of  subsets  of  X,  where  Mi  =  C(h^),  the  goal  is  to  produce  a 
sequence 

n2  ^  ...  ^  m  =  y 

of  subsets  of  X  that  somehow  summarize  or  represent  the 
underlying  picture  that  led  to  the  original  input  data. 

In  that  depends  only  upon  this  amounts  to  defining  an 

isotone  mapping  y  on  <?(X)  . 

Not  every  isotone  mapping  on  (?(X)  is  appropriate  for 
defining  a  flat  segmentation  method.  One  wants  to  attach 
some  spatial  significance  to  the  decision  as  to  whether  a 
point  x  belongs  to  y (M) .  One  way  of  doing  this  is  to  insist 
that  y  be  point-based  in  that  for  each  x  e  y(X)  there  is  a 

subset  N(x)  of  X  containing  x  such  that: 

(PBl )  y(0)  =  0,  and  y(N(x))  ?  0. 

(PB2)  For  A£  N(x),  if  y  (A)  ¥■  0 ,  then  x  e  y  (A)  . 

(PB3)  For  M  C  x,  x  e  y (M)  if  and  only  if  x  e  y (MON (x) ) . 

It  is  immediate  that  for  every  subset  M  of  X, 

Y(M)  =  (J  y(MhN(x}).  If  one  wants  the  output  of  y  to 
x  e  y (X) 

be  independent  of  the  polarity  of  the  image,  it  is  useful  to 
also  insist  that  Y  preserve  complements  in  that  it  also  satisfy 

(PB4)  For  M  c  X,  y  (X\M)  =  Y(X)\Y(M). 

The  local  version  of  this  is  the  content  of  Theorem  1. 


I 
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Theorem  1.  For  a  point-based  isotone  mapping  V  on  G>(x) , 
axiom  (PB4 )  is  equivalent  to  the  assertion  that  for  each 
subset  A  of  N(x)  ,  exactly  one  of  Y (A)  and  y(N(xNA)  shall 
be  nonempty. 

! 

If  (n(x)}  and  {M(x) }  are  each  families  of  neighborhoods 
satisfying  axioms  (PB1)  through  (PB3)  ,  then  so  is  {M(x)f\N(x)}. 
In  view  of  this  there  is  no  harm  in  assuming; 

I 

(P35)  If  M (x)  satisfies  (PBl)  through  (PB3) ,  then 
N (x)  C.  M{x)  . 


Such  a  minimal  family  of  subsets  of  X  will  be  called  tne 
system  of  neighborhoods  associated  with  y.  Henceforth,  when 


we  speak  of  a  point-based  isotone  mapping  Y,  it  will  always 
be  assumed  that  {N(x)>x  £  y  (x)  denotes  the  system  of  neigh¬ 
borhoods  associated  with  y. 


Definition.  The  point-based  isotone  mapping  y  is 
said  to  be  frequency-defined  if  ther  exist  positive  integers 

j  and  k  such  that  for  every  x  e  y  (X)  ,  (i)  k  =  #W(x),  and 
(ii)  x  e  y(M)  if  and  only  if  j  <  #(Mr»N(x)).  Here  # A  denotes 
the  cardinality  of  the  set  A. 


Theorem  2.  Let  y  be  a  point-based,  complement-pre¬ 
serving  isotone  mapping  on  CKX) .  Necessary  and  sufficient 
conditions  for  y  to  be  frequency-defined  are  that; 

(i)  the  neighborhoods  N(x)  all  have  the  same  cardinality 
k,  where  k  is  odd,  and 

(ii)  if  has  minimal  cardinality  among  those  subsets 


i 


A  of  N (x)  for  which  x  e  y  (A)  ,  then 


#AX  =  (1  +  k)/2. 

To  say  that  the  mapping  y  is  a  join  homomorphism  is  to 
say  that  y (MUN)  =  y (M)  U  y (N) .  There  is  a  dual  notion  of 
meet  homomorphism  and  to  say  that  y  is  a  homomorphism  is  to 

say  that  it  is  both  a  join  and  a  meet  homomorphism.  The 
next  theorem  shows  that  these  conditions  are  rather 
powerful,  and  will  only  be  met  in  trivial  situations. 

Theorem  3.  Let  y  be  a  point-based  isotone  mapping 
on  (?(X)  .  Then  (la)  «=»  (lb)  ,  (2a)  <=>  (2b)  and  (3aK=»  (3b)  . 

(la)  y  is  a  join  homomorphism. 

(lb)  For  each  x  e  y(X)  there  is  a  subset  A.,  of  N(x) 
such  that  x  e  y  (M)  if  and  only  if  Mr>A..  /  0. 

(2a)  y  is  a  meet  homomorphism. 

(2b)  Corresponding  to  each  x  r.  y(X)  there  is  a  subset 
N(x)  such  that  x  e  y  (M)  if  and  only  if  Q  M. 

(3a)  y  is  a  homomorphism. 

(3b)  For  each  x  e  y(X)  there  corresponds  an  element 

yx  of  N(x)  such  that  x  s  y  (M)  if  and  only  if  y^_  e  M. 

Corollary  4.  If  y  is  also  complement-preserving, 
then  all  six  conditions  of  the  theorem  are  mutually  equivalent 
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A  translation  on  X  is  a  mapping  of  the  form  where 

p,q  are  fixed  integers  and  T  (i,j)  =  (p+i,q+j).  Unless 

P/Q 

p  =  0  =  q,  the  domain  and  image  of  T  will  be  proper  subsets 

p ,  q 

of  X.  We  agree  to  call  the  point-based  isotone  mapping 
translation- invariant  provided  it  satisfies: 

(P36)  If  N(x)  is  contained  in  the  domain  of  the 

translation  T,  and  if  y  =  T(x),  then  N(y)  = 

T  (N  (x)  )  ,  and  for  AC  N(x),  x  e  (A)  if  and 
only  if  y  e  y (T (A) ) . 

Remark  5.  If  M  is  acted  upon  by  a  translation  to 
produce  N,  it  is  desirable  that  y (M)  produce  y (N)  when  it  is 
acted  on  by  that  same  translation.  As  long  as  y  is  translation- 
invariant,  MS  domain  T,  and  MU  T(M)  C  y  (X)  ,  it  is  easy  to 
see  tnat  this  is  true.  An  example  is  provided  in  [8J  to 
show  that  this  can  fail  even  if  y  is  translation-invariant. 

Remark  6.  The  actual  construction  of  a  point-based 
isotone  mapping  on  (?(X)  may  now  easily  be  understood.  A 
system  of  neighborhoods  {N(x)}  for  points  x  in  some  subset  Y 
of  X  is  first  chosen.  One  then  defines  a  family  of  mappings 
(yx) x  E  y'  w*iere  Yx  ^  isotone  mapping  on  N(x)  such  that: 

(1)  y  (0)  =  0,  and  y  (N (x) )  =  {x}. 

A  X 

To  produce  a  complement-preserving  mapping,  one  also  wants 

(2)  For  AC  n(x)  ,  exactly  one  of  y^  (A)  and  yx(N(x)^A) 
is  nonempty. 

The  mapping  y  is  now  defined  by  the  rule  x  e  y (M)  if  and  only 
if  yx(MAN(x))  =  {x}.  This  is  the  technique  that  will  be 

used  for  the  remainder  of  the  paper. 
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IV.  Underlying  statistical  considerations. 

The  construction  of  a  flat  segmentation  method  involves 

the  definition  of  an  isotone  mapping  y  on  (r>(X) .  The  input 

data  represents  a  subset  M  of  X  that  is  an  estimate  of  the 
* 

subset  M  of  X  that  is  the  true  input  data.  The  idea  is 

to  try  to  define  y  so  that  y (M)  is  in  some  sense  a  better 

* 

estimate  of  M  than  is  M.  A  crude  statistical  model  may  be 
constructed  by  assuming  that  membership  in  M  has  probability 

p  (p  >  0.5)  of  providing  a  correct  estimate  of  membership  in 

*  ★ 

M  ,  that  membership  in  XvM  has  a  probability  q  (q  >  0.5) 

* 

of  providing  a  correct  estimate  of  membership  in  XvM  ,  and 

that  these  probabilities  are  independently  distributed  over 

* 

the  members  of  X.  An  interior  point  of  14  is  defined  to  be 

* 

a  point  x  for  which  N(x)  C  m  ;  dually,  x  is  an  exterior 
*  * 

point  of  M  if  N(x) £  X^M  .  The  gain  of  the  isotone  mapping 

y  is  defined  to  be  the  sum  of  the  probability  of  correct 

★ 

classification  for  an  interior  point  of  M  and  that  of  an 

exterior  point.  Needless  to  say,  these  are  idealized  concepts, 

* 

as  M  is  precisely  the  unknown  set  that  we  are  trying  to 
estimate.  If  G(y)  denotes  the  gain  f  y,  we  tnen  have 

Lemma  7.  Let  y  be  complement-preserving  with  k  =  #N(x) 


system  N(x)  of  neighborhoods  with  y'  frequency-defined,  then 


G(y)  *  G(y  ' ) 
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When  p  =  q,  the  above  theorem  shows  that  the  "gain" 

from  a  flat  segmentation  method  can  be  maximized  by  using  a 

point-based  isotone  mapping  that  is  complement-preserving 

and  frequency-defined.  With  this  thought  in  mind,  we  shall 

concentrate  on  such  mappings,  defining  them  using  k  by  k 

2 

neighborhoods  of  points  with  k  an  odd  integer.  The  j/k 
rule  will  be  the  unique  such  mapping  defined  on  a  k  by  k 
neighborhood,  and  the  3/5  rule  will  refer  to  the  y  that  is 
based  on  a  point  together  with  its  4  immediate  direct  neighbors 
(North,  Soucn,  East  and  West).  These  rules  were  discussed 
in  some  detail  in  f7]  and  [8] ,  and  examples  were  given  there 

•k 

to  show  that  near  the  boundary  of  M  ,  the  use  of  these  rules 
can  actually  decrease  the  probability  of  correct  classifi¬ 
cation.  This  suggests  using  (see  0T]  ,  Table  3)  an  isotone 
mapping  y  that  is  based  on  a  weighted  mean  or  on  M^N(x) 
containing  certain  desirable  subsets.  For  example,  one 
would  be  more  likely  to  conclude  that  x  e  y (A)  if  A  were 
the  left-hand  subset  than  if  A  were  the  right-hand  set: 


0  0  0  0  0 
0  1110 
0  1110 
0  1110 
0  0  0  0  0 


10  0  10 
0  0  0  0  1 
10  0  10 
0  10  0  1 
1  0  0  0  1 


V.  The  SEGMENT  Program. 

This  was  described  in  some  detail  in  C63 ,  T7]  and  to  a 
lesser  extent  in  [8j .  For  the  reader's  convenience,  a  brief 
description  is  also  included  here. 


L  3' 


* 

! 

1 


h 


i 


l 


rl 


A.  Input .  An  m  by  n  matrix  A  having  nonnegative 
integers  as  entries. 

B.  Pre filtering.  A  k  by  k  mean  or  median  filter  is 

applied  to  smooth  the  data  of  A,  and  the  output  is  rounded 
to  produce  B.  Spot  noise  is  removed  by  deleting  the  highest 
and  lowest  i  values  when  computing  the  k  by  k  mean. 

C.  Tiiinout.  A  frequency  count  is  made  of  the  values 

appearing  in  B.  Those  values  that  occur  with  frequency  less 
than  some  thresnhold  are  deleted,  and  a  nearest  value  rule 
used  to  reassign  them.  Alternately,  the  k  highest  occuring 
values  can  be  retained.  The  resulting  matrix  is  denoted  C. 

D.  Suppose  matrix  C  has  data  values  j ^ ,  j2,  ...,  j^. 

The  nondirected  3  by  3  dispersion  of  value  j1  is  defined  as 

follows:  For  each  point  in  C  having  value  j ^ ,  look  at  a 

3  by  3  neignborhood  centered  on  that  point.  Calculate  the 
number  of  points  in  that  neighborhood  that  do  not  have  value 
j^,  and  average  this  figure  over  alJ  points  having  value  j^. 

This  is  the  dispersion  for  cluster  j^.  One  can  do  a  similar 
calculation  for  j2,  . ..,  j^.  An  analogous  definition  can  be 

formulated  for  larger  neighborhoods,  and  certain  options  are 
included  that  compute  directional  versions  of  the  dispersion. 

E.  Various  options  exist  for  deleting  one  or  more  of 
the  data  values  having  high  dispersion  levels.  The  simplest 
is  to  just  delete  the  highest  value  dispersion.  Other 
options  might  include  deleting  all  data  values  wnose  dis¬ 
persion  exceeds  0.9  on  the  first  pass,  and  then  dropping 
this  cutoff  down  by  increments  of  0.1  until  a  stopping 
criterion  is  reacned.  Various  choices  also  exist  for  the 
reassignment  of  points  whose  value  has  been  deleted.  They 
proceed  on  both  a  global  and  a  local  basis.  Cluster  means 


.W 
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can  be  computed  for  each  remaining  cluster,  and  pixels 
reassigned  to  tne  cluster  to  which  their  input  value  is 
closest,  or  an  entire  cluster  can  be  reassigned  to  the 
cluster  to  which  its  mean  is  closest.  A  second  technique 
involves  the  reassignment  of  points  to  the  next  higher  or 
lower  cluster  that  is  still  valid,  either  on  a  local  or  a 
global  basis.  As  these  techniques  have  been  described, 
they  are  not  monotone  equivariant,  but  if  the  data  are 
first  rank  ordered,  and  the  output  is  labeled  according  to 
which  value  corresponds  to  which  rank,  the  resulting 
techniques  do  become  monotone  equivariant. 


VI  The  SLICER  program. 

A  crude  segmentation  program  can  easily  be  devised 
using  the  material  of  Sections  III  and  IV.  It  has  two 
phases  s 

Phase  1.  Find  levels  at  which  to  slice  the  data. 

2 

Phase  2 .  Make  appropriate  si '  .es  and  use  a  j/k  rule 
to  increase  the  probability  of  correct  classification. 

The  Phase  2  portion  of  the  program  is  a  flat  segmentation 
method.  The  program  itself  can  either  be  used  as  a  cleaning 
algorithm,  or  as  a  tool  for  obtaining  a  first  approximation 
to  the  principal  regions  of  a  digital  picture. 

Description  of  Phase  1.  A  total  of  six  methods  of 
determining  the  slice  levels  are  built  into  the  program. 

Method  1 .  Compute  the  mean  M  and  standard  deviation 

S  of  the  input  data,  and  take  slices  at  M  -  2S,  M  -  S,  M, 

M  +  S,  and  M  +  2S .  If  the  extreme  values  are  outside  the 
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data  range,  they  are  modified  to  bring  them  within  the 
minimum  and  maximum  values  of  the  data. 

Method  2 .  Slice  at  relative  maxima  of  the  histogram 
of  the  data. 

Method  3 .  Slice  at  relative  minima  of  the  nondirected 
3  by  3  dispersion. 

Method  4 .  This  is  similar  to  Method  3,  except  that 
if  A  is  the  vector  of  distinct  data  values,  then  at  level 
A[i]  ,  one  is  interested  in  those  values  that  lie  outside 
the  range  from  A[i  -  1]  to  A[i  +  1]  . 

Method  5 .  Slice  at  relative  minima  of  the  average 
variance  in  t  by  t  neighborhoods. 

Method  6.  User  inputs  slice  values. 

description  of  Phase  2.  Having  arrived  at  a  list  of 
2 

slice  levels,  a  j/k  rule  is  chosen.  The  actual  slices  are 
constructed  so  that  they  are  halfway  Detween  the  slice  input 
levels.  Thus  if  one  wanted  slices  at  17,  20,  22,26,  and  30, 
then  one  would  look  at  the  sets  M^  c  m2  c  M3  <L  Il4  C  Mg, 

where  fK  =  {x:  x  <  k^}  with  k^  =  18.5,  k2  =  21,  k3  =  24, 

k.  =  28  and  kc  =  the  maximum  value  of  X.  This  forms  5 
4  5 

clusters,  and  they  are  assigned  the  values  of  17,  20,  22, 

26  and  30.  There  are  also  occasions  when  one  would  want  the 

actual  slice  values  to  coincide  with  the  input  slice  values. 

2 

If  y  represents  the  isotone  mapping  specified  by  the  j/k 
rule,  the  output  would  then  be  the  sequence  ; 


y(M1)C  Y(fi2)  c  y<m3)  c.  y(m4)  C  y(M5). 


If  the  actual  slices  are  taken  at  levels  k2,...,kt 
the  implementation  of  this  process  is  very  simple.  One 
starts  by  forming  a  Boolean  matrix  with  1  denoting  a  value 
<  k^  and  0  otherwise.  One  then  looks  at  k  by  k  regions. 

If  there  are  j  or  more  l’s  in  such  a  region,  the  output 
is  1;  otherwise  it  is  0.  The  output  for  level  k2  is  added 
to  this,  and  the  process  continues  until  the  supply  of 
slice  levels  is  exhausted. 

The  actual  implementation  of  the  SLICER  program  can 

be  done  more  efficiently  by  first  performing  a  k  by  k 

median  filter,  and  then  applying  the  algorithm  described 

at  the  beginning  of  Pnase  2.  That  way,  if  one  wanted  to 

examine  the  effect  of  various  types  of  slicing,  one  would 

2 

only  have  to  perform  the  j/k  rule  a  single  time,  as  this 
decision  rule  is  in  reality  implemented  by  the  median  filter. 


VI I  Some  Examples. 

Plates  1-4  contain  a  comparison  of  the  output  of  the 
SEGMENT  program  with  that  of  prograir  SLICER.  The  plates 
all  follow  the  same  format,  so  a  single  explanation  of 
them  should  suffice.  There  are  6  pictures  per  plate.  The 
upper  left  picture  is  the  original  data  set.  The  middle 
picture  on  the  upper  level  is  the  SEGMENT  output  using  a 
nondirected  3  by  3  dispersion,  a  3  by  3  mean  filter,  and  the 
version  of  the  program  that  removes  1  cluster  at  a  time, 
reassigning  its  members  locally  to  the  next  higher  or  lower 
available  cluster.  The  parameters  were  set  to  produce  5 
clusters.  The  remaining  pictures  relate  to  the  SLICER 
program,  though  it  should  be  noted  that  the  program  was 
modified  to  make  it  iterate  with  a  stopping  criterion  set 
for  when  an  iteration  did  not  change  the  cluster  means. 

At  each  stage,  cluster  means  were  computed,  and  these  values 
used  as  slice  level  input  for  the  next  iteration.  The  upper 
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right  picture  was  the  result  of  the  first,  third  and  fifth 
levels  of  option  1  of  SLICER,  the  bottom  left  was  the  result 
of  option  1,  the  bottom  middle  the  result  of  3  slices 
equally  placed  between  the  minimum  and  maximum  values  of 
tne  input  data,  and  the  bottom  right  the  result  of  5  slice 
levels  equally  placed.  Thus  for  3  levels,  the  slices  are 
taken  1/4,  2/4  and  3/4  of  the  way  from  the  minimum  to  the 
maximum  data  value;  for  5  levels,  the  slices  are  1/6,  2/6, 

3/6,  4/6  and  5/6  of  the  way.  In  each  case  a  5/9  rule  was 
applied  to  produce  the  final  output.  No  further  cleaning 
algorithm  was  applied.  It  should  be  noted  that  the  various 
SLICER  outputs  produced  reasonable  crude  segmentations  of 
the  data.  Further  work  is  clearly  needed  on  the  initial 
selection  of  slice  levels.  The  outputs  are  especially 
useful  when  one  resorts  to  line  printer  representations, 
as  one  can  intelligently  reduce  a  picture  to  a  desired 
number  of  grey  levels. 

Plate  3  is  worthy  of  note.  It  illustrates  one  of  the 
pitfalls  of  the  SEGMENT  program.  This  program  operates 
essentially  by  merging  clusters;  no  provision  is  made  for 
a  subsequent  splitting  of  a  cluster,  :nless  the  dispersion 
criterion  deems  that  cluster  to  be  too  noisy  to  be  maintained. 
The  vehicle  in  the  center  of  the  picture  is  too  large 
because  of  some  early  mergers  that  were  not  undone.  In 
the  lower  rignt-hand  corner  of  Plate  4,  a  portion  of  the 
Gulf  Stream  is  visible.  In  all  but  1  of  the  outputs,  this 
cluster  was  merged  with  warm  coastal  water.  This  was 
caused  largely  by  the  iterations  using  cluster  means  as 
outputs.  The  original  iteration  more  correctly  merged  the 
Gulf  Stream  into  the  same  cluster  as  the  land  masses, 
tnus  correctly  identifying  it  as  being  much  warmer  than 
the  coastal  water. 


Portion  of  Record  12  of  Westingnouse  FLIR 
data  tape. 


Portion  of  weather  tape  furnisned  by  Air  Force 
Geophysical  Laboratory.  The  picture  shows  the 
region  off  of  the  Chesapeake  Bay  with  the  edge 
of  the  Gulf  Stream  visible  in  the  lower  right 


corner. 
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